Load packages

library(tidyverse)
library(gridExtra)
library(cowplot)
library(ggridges)
library(ggstance)
library(treeio)
library(ggtree)
library(tidytree)

Load and combine data

species <- read_csv("../data/species_data.csv") %>% select(-species)
data <- read_csv("../data/data.csv") %>% 
  left_join(species, by = c("species" = "species_formatted")) %>% 
  rename(label = species_latin)
data2 <- data %>% group_by(species, label) %>% summarise(n = sum(n))
fulltree <- read.nexus("../data/consensusTree_10kTrees_298Primates_V3.nex")
refs <- read_csv("../data/ref_nodes.csv")
# turn tree into tidy dataframe
tree2 <- as_tibble(fulltree)

tree3 <- tree2 %>% 
  left_join(data2) %>% 
  mutate(
    hasN = ifelse(is.na(n), 0, .5),
    hasN2 = ifelse(is.na(n), .1, .5)) %>% 
  left_join(refs) %>% 
  # # also merge w datasheet listing node, n for inner nodes
  # left_join(Ns, by = "node") %>% 
  # mutate(n = coalesce(n.x, n.y)) %>% 
  # select(-n.x, -n.y) %>% 
  groupClade(c(493, 496, 429, 302, 408)) %>% 
  mutate(group = fct_recode(group, "2" = "1"))
Joining, by = "label"
Joining, by = "node"
# turn back into tree
tree4 <- as.treedata(tree3)

Figure out nodes

This makes a rectangular and a circular tree with the node numbers displayed for reference (saved in the graphs folder).

tree3.2 <- as.treedata(tree3)
# display node numbers for reference
ggtree(tree3.2) +
  # tip labels
  geom_tippoint(aes(size = n), col = "seagreen", alpha = .5) +
  geom_tiplab(offset = 1, size = 3) +
  # node labels
  geom_text(aes(label = node, x = branch), size = 2, col = "blue", vjust = -.5) +
  expand_limits(x = 90) +
  # display timescale at the bottom
  theme_tree2()
ggsave("../graphs/full_tree_nodes.pdf", width = 8, height = 20, scale = 2)
ggtree(tree3.2, layout = "circular") +
  geom_tippoint(aes(size = n), col = "seagreen", alpha = .5) +
  geom_tiplab2(offset = 2, size = 3) +
  geom_text2(aes(label = node), size = 1.5, col = "blue") +
  xlim(NA, 100)
ggsave("../graphs/full_tree_nodes_circular.pdf", width = 8, height = 8, scale = 2)

Circular tree of 298 primates

cols <- viridis::viridis(4, end = .9)
# base plot
p <- ggtree(tree4, aes(size = hasN, alpha = hasN2), layout = "circular") +
  # root
  geom_rootpoint(size = 1) +
  # tips
  geom_tippoint(aes(size = n), alpha = .5) +
  geom_tiplab2(aes(alpha = hasN), offset = 2, size = 3) +
  # tweak scales
  scale_alpha_continuous(range = c(.3, 1)) +
  scale_size_continuous(range = c(.5, 8)) +
  # widen plotting area
  xlim(NA, 100)
# highlight clades with background colors
p2 <- p + 
  geom_hilight(node = 493, fill = cols[1], alpha = .3) +
  geom_hilight(node = 496, fill = cols[1], alpha = .3) +
  geom_hilight(node = 429, fill = cols[2], alpha = .3) +
  geom_hilight(node = 303, fill = cols[3], alpha = .3) +
  geom_hilight(node = 408, fill = cols[4], alpha = .3) +
  # plot tree again to be on top of the highlights
  geom_tree() +
  geom_rootpoint(size = 1)
# highlight clades with branch colors
p3 <- p + 
  aes(col = group) +
  scale_color_manual(values = c("gray30", cols))
plots <- mget(c("p", "p2", "p3"))
grid.arrange(p, p2, p3, nrow = 1)

# png with 3x1
ggsave("../graphs/phylo_full.png", arrangeGrob(grobs = plots, nrow = 1), 
       width = 24, height = 8, scale = 2, dpi = 72)
Removed 287 rows containing missing values (geom_point_g_gtree).Removed 287 rows containing missing values (geom_point_g_gtree).Removed 287 rows containing missing values (geom_point_g_gtree).
# pdf with 1 per page
ggsave("../graphs/phylo_full.pdf", marrangeGrob(grobs = plots, nrow = 1, ncol = 1), 
       width = 8, height = 8, scale = 2, dpi = 72)
Removed 287 rows containing missing values (geom_point_g_gtree).Removed 287 rows containing missing values (geom_point_g_gtree).Removed 287 rows containing missing values (geom_point_g_gtree).

Sample size in detail

# subset tree to just those species who have sample sizes reported, i.e. those who were tested
to_drop <- tree3 %>% filter(is.na(n)) %>% pull(label)
tree5 <- drop.tip(tree4, to_drop)
data3 <- data %>% 
  select(label, everything()) %>% 
  rename(num = n)

# species with more than X sites can get a density?
d3a <- data3 %>% group_by(species) %>% filter(n_distinct(site) >= 4)
d3b <- data3 %>% # setdiff(data3, d3a) ## <- do setdiff instead to NOT show points for densities
  group_by(species) %>% 
  # create variable num2 is NA if there's only one data point for a species
  # --> those species will only get the vertical crossbar
  mutate(flag = n_distinct(site) == 1) %>% 
  ungroup %>% 
  mutate(num2 = ifelse(flag, NA, num))

# for vertical crossbar = median
data4 <- data3 %>% 
  group_by(label, species) %>% 
  summarise(Mdn = median(num)) # totalN = sum(num), sitesN = n_distinct(site)
q <- ggtree(tree5, aes(col = group)) +
  # this is a dummy point to expand the x scale
  # the typical ways weirdly also expand it for the sample size panel once that's added
  geom_point(data = tibble(x = 135, y = 1, .panel = "Tree", group = NA)) +
  # tip labels
  geom_tippoint(aes(size = n), alpha = .5) +
  geom_tiplab(offset = 4, size = 3) +
  # tweak scales
  scale_color_manual(values = c("grey30", cols)) +
  scale_fill_manual(values = cols[4]) + # when all categories are taken: cols
  # display timescale at the bottom
  theme_tree2()
# trying out which geoms to layer
s1 <- ggplot(data3, aes(x = num, y = label, group = label)) +
  geom_density_ridges(lwd = .3, col = "grey80", bandwidth = 1)

s2 <- ggplot(data3, aes(x = num, y = label, group = label)) +
  geom_crossbarh(data = data4, aes(x = Mdn, xmin = Mdn, xmax = Mdn), width = .5) +
  geom_point(alpha = .5)

s3 <- ggplot(data3, aes(x = num, y = label, group = label)) +
  geom_boxploth()

plot_grid(s1, s2, s3, nrow = 1)

# right-side viz depends on the number of sites per species:
# 1 site = vertical crossbar only
# 2+ sites = points + crossbar at median
# X+ sites = densities (currently, X = 4 just to illustrate)

# dirty hack: there's a small character in front of "Sample sizes" to have that panel sort to the right (alphabetically) until I figure out why it doesn't just go by order. This cropped up as an issue when I added the dummy point for the x-axis expansion...

q %>% 
  facet_plot("Ù´ Sample sizes", 
             d3a, geom_density_ridges, aes(x = num, group = label, fill = group), 
             alpha = .5, lwd = .3) %>% 
  facet_plot("Ù´ Sample sizes",
             data4, geom_crossbarh, aes(x = Mdn, xmin = Mdn, xmax = Mdn, group = label, 
             col = group), alpha = .5, width = .3, position = position_nudge(y = -.1)) %>%
  facet_plot("Ù´ Sample sizes", 
             d3b, geom_point, aes(x = num2, group = label), alpha = .5, 
             position = position_nudge(y = -.1)) +
  panel_border()

ggsave("../graphs/phylo_ridge.png", width = 4, height = 3, scale = 2)

Relevant functions to tinker with, if necessary:

facet_plot
function (p, panel, data, geom, mapping = NULL, ...) 
{
    p <- add_panel(p, panel)
    df <- p %+>% data
    p + geom(data = df, mapping = mapping, ...)
}
<bytecode: 0x7f8a0db5f978>
<environment: namespace:ggtree>
ggtree:::add_panel
function (p, panel) 
{
    df <- p$data
    if (is.null(df[[".panel"]])) {
        df[[".panel"]] <- factor("Tree")
    }
    levels(df$.panel) %<>% c(., panel)
    p$data <- df
    p + facet_grid(. ~ .panel, scales = "free_x")
}
<bytecode: 0x7f8a0db62660>
<environment: namespace:ggtree>

Session info

sessionInfo()
R version 3.6.1 (2019-07-05)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS Mojave 10.14.5

Matrix products: default
BLAS:   /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] tidytree_0.2.8  ggtree_1.16.6   treeio_1.8.2    ggstance_0.3.3  ggridges_0.5.1 
 [6] cowplot_1.0.0   gridExtra_2.3   forcats_0.4.0   stringr_1.4.0   dplyr_0.8.3    
[11] purrr_0.3.2     readr_1.3.1     tidyr_1.0.0     tibble_2.1.3    ggplot2_3.2.1  
[16] tidyverse_1.2.1

loaded via a namespace (and not attached):
 [1] tidyselect_0.2.5   xfun_0.10          reshape2_1.4.3     haven_2.1.1        lattice_0.20-38   
 [6] colorspace_1.4-1   vctrs_0.2.0        generics_0.0.2     viridisLite_0.3.0  rlang_0.4.0       
[11] pillar_1.4.2       glue_1.3.1         withr_2.1.2        modelr_0.1.5       readxl_1.3.1      
[16] rvcheck_0.1.5      lifecycle_0.1.0    plyr_1.8.4         munsell_0.5.0      gtable_0.3.0      
[21] cellranger_1.1.0   rvest_0.3.4        labeling_0.3       knitr_1.25         parallel_3.6.1    
[26] broom_0.5.2        Rcpp_1.0.2         BiocManager_1.30.7 scales_1.0.0       backports_1.1.5   
[31] jsonlite_1.6       hms_0.5.1          stringi_1.4.3      grid_3.6.1         cli_1.1.0         
[36] tools_3.6.1        magrittr_1.5       lazyeval_0.2.2     crayon_1.3.4       ape_5.3           
[41] pkgconfig_2.0.3    zeallot_0.1.0      xml2_1.2.2         lubridate_1.7.4    viridis_0.5.1     
[46] assertthat_0.2.1   httr_1.4.1         rstudioapi_0.10    R6_2.4.0           nlme_3.1-141      
[51] compiler_3.6.1    
LS0tCnRpdGxlOiAiUGh5bG9nZW5ldGljIFRyZWUiCm91dHB1dDogCiAgaHRtbF9ub3RlYm9vazoKICAgIGNzczogc3R5bGUuY3NzCiAgICB0aGVtZTogcGFwZXIKLS0tCgpMb2FkIHBhY2thZ2VzCgpgYGB7ciwgbWVzc2FnZT1GQUxTRX0KbGlicmFyeSh0aWR5dmVyc2UpCmxpYnJhcnkoZ3JpZEV4dHJhKQpsaWJyYXJ5KGNvd3Bsb3QpCmxpYnJhcnkoZ2dyaWRnZXMpCmxpYnJhcnkoZ2dzdGFuY2UpCmxpYnJhcnkodHJlZWlvKQpsaWJyYXJ5KGdndHJlZSkKbGlicmFyeSh0aWR5dHJlZSkKYGBgCgpMb2FkIGFuZCBjb21iaW5lIGRhdGEKCmBgYHtyLCBtZXNzYWdlPUZBTFNFfQpzcGVjaWVzIDwtIHJlYWRfY3N2KCIuLi9kYXRhL3NwZWNpZXNfZGF0YS5jc3YiKSAlPiUgc2VsZWN0KC1zcGVjaWVzKQpkYXRhIDwtIHJlYWRfY3N2KCIuLi9kYXRhL2RhdGEuY3N2IikgJT4lIAogIGxlZnRfam9pbihzcGVjaWVzLCBieSA9IGMoInNwZWNpZXMiID0gInNwZWNpZXNfZm9ybWF0dGVkIikpICU+JSAKICByZW5hbWUobGFiZWwgPSBzcGVjaWVzX2xhdGluKQpkYXRhMiA8LSBkYXRhICU+JSBncm91cF9ieShzcGVjaWVzLCBsYWJlbCkgJT4lIHN1bW1hcmlzZShuID0gc3VtKG4pKQpmdWxsdHJlZSA8LSByZWFkLm5leHVzKCIuLi9kYXRhL2NvbnNlbnN1c1RyZWVfMTBrVHJlZXNfMjk4UHJpbWF0ZXNfVjMubmV4IikKcmVmcyA8LSByZWFkX2NzdigiLi4vZGF0YS9yZWZfbm9kZXMuY3N2IikKYGBgCgpgYGB7cn0KIyB0dXJuIHRyZWUgaW50byB0aWR5IGRhdGFmcmFtZQp0cmVlMiA8LSBhc190aWJibGUoZnVsbHRyZWUpCgp0cmVlMyA8LSB0cmVlMiAlPiUgCiAgbGVmdF9qb2luKGRhdGEyKSAlPiUgCiAgbXV0YXRlKAogICAgaGFzTiA9IGlmZWxzZShpcy5uYShuKSwgMCwgLjUpLAogICAgaGFzTjIgPSBpZmVsc2UoaXMubmEobiksIC4xLCAuNSkpICU+JSAKICBsZWZ0X2pvaW4ocmVmcykgJT4lIAogICMgIyBhbHNvIG1lcmdlIHcgZGF0YXNoZWV0IGxpc3Rpbmcgbm9kZSwgbiBmb3IgaW5uZXIgbm9kZXMKICAjIGxlZnRfam9pbihOcywgYnkgPSAibm9kZSIpICU+JSAKICAjIG11dGF0ZShuID0gY29hbGVzY2Uobi54LCBuLnkpKSAlPiUgCiAgIyBzZWxlY3QoLW4ueCwgLW4ueSkgJT4lIAogIGdyb3VwQ2xhZGUoYyg0OTMsIDQ5NiwgNDI5LCAzMDIsIDQwOCkpICU+JSAKICBtdXRhdGUoZ3JvdXAgPSBmY3RfcmVjb2RlKGdyb3VwLCAiMiIgPSAiMSIpKQoKIyB0dXJuIGJhY2sgaW50byB0cmVlCnRyZWU0IDwtIGFzLnRyZWVkYXRhKHRyZWUzKQpgYGAKCiMgRmlndXJlIG91dCBub2RlcwoKVGhpcyBtYWtlcyBhIHJlY3Rhbmd1bGFyIGFuZCBhIGNpcmN1bGFyIHRyZWUgd2l0aCB0aGUgbm9kZSBudW1iZXJzIGRpc3BsYXllZCBmb3IgcmVmZXJlbmNlIChzYXZlZCBpbiB0aGUgYGdyYXBoc2AgZm9sZGVyKS4KCmBgYHtyfQp0cmVlMy4yIDwtIGFzLnRyZWVkYXRhKHRyZWUzKQpgYGAKCmBgYHtyLCBmaWcud2lkdGg9OCwgZmlnLmhlaWdodD0yMCwgZXZhbD1GQUxTRX0KIyBkaXNwbGF5IG5vZGUgbnVtYmVycyBmb3IgcmVmZXJlbmNlCmdndHJlZSh0cmVlMy4yKSArCiAgIyB0aXAgbGFiZWxzCiAgZ2VvbV90aXBwb2ludChhZXMoc2l6ZSA9IG4pLCBjb2wgPSAic2VhZ3JlZW4iLCBhbHBoYSA9IC41KSArCiAgZ2VvbV90aXBsYWIob2Zmc2V0ID0gMSwgc2l6ZSA9IDMpICsKICAjIG5vZGUgbGFiZWxzCiAgZ2VvbV90ZXh0KGFlcyhsYWJlbCA9IG5vZGUsIHggPSBicmFuY2gpLCBzaXplID0gMiwgY29sID0gImJsdWUiLCB2anVzdCA9IC0uNSkgKwogIGV4cGFuZF9saW1pdHMoeCA9IDkwKSArCiAgIyBkaXNwbGF5IHRpbWVzY2FsZSBhdCB0aGUgYm90dG9tCiAgdGhlbWVfdHJlZTIoKQpgYGAKCmBgYHtyLCBldmFsPUZBTFNFfQpnZ3NhdmUoIi4uL2dyYXBocy9mdWxsX3RyZWVfbm9kZXMucGRmIiwgd2lkdGggPSA4LCBoZWlnaHQgPSAyMCwgc2NhbGUgPSAyKQpgYGAKCmBgYHtyLCBmaWcud2lkdGg9OCwgZmlnLmhlaWdodD04LCBldmFsPUZBTFNFfQpnZ3RyZWUodHJlZTMuMiwgbGF5b3V0ID0gImNpcmN1bGFyIikgKwogIGdlb21fdGlwcG9pbnQoYWVzKHNpemUgPSBuKSwgY29sID0gInNlYWdyZWVuIiwgYWxwaGEgPSAuNSkgKwogIGdlb21fdGlwbGFiMihvZmZzZXQgPSAyLCBzaXplID0gMykgKwogIGdlb21fdGV4dDIoYWVzKGxhYmVsID0gbm9kZSksIHNpemUgPSAxLjUsIGNvbCA9ICJibHVlIikgKwogIHhsaW0oTkEsIDEwMCkKYGBgCgpgYGB7ciwgZXZhbD1GQUxTRX0KZ2dzYXZlKCIuLi9ncmFwaHMvZnVsbF90cmVlX25vZGVzX2NpcmN1bGFyLnBkZiIsIHdpZHRoID0gOCwgaGVpZ2h0ID0gOCwgc2NhbGUgPSAyKQpgYGAKCiMgQ2lyY3VsYXIgdHJlZSBvZiAyOTggcHJpbWF0ZXMKCmBgYHtyfQpjb2xzIDwtIHZpcmlkaXM6OnZpcmlkaXMoNCwgZW5kID0gLjkpCmBgYAoKYGBge3J9CiMgYmFzZSBwbG90CnAgPC0gZ2d0cmVlKHRyZWU0LCBhZXMoc2l6ZSA9IGhhc04sIGFscGhhID0gaGFzTjIpLCBsYXlvdXQgPSAiY2lyY3VsYXIiKSArCiAgIyByb290CiAgZ2VvbV9yb290cG9pbnQoc2l6ZSA9IDEpICsKICAjIHRpcHMKICBnZW9tX3RpcHBvaW50KGFlcyhzaXplID0gbiksIGFscGhhID0gLjUpICsKICBnZW9tX3RpcGxhYjIoYWVzKGFscGhhID0gaGFzTiksIG9mZnNldCA9IDIsIHNpemUgPSAzKSArCiAgIyB0d2VhayBzY2FsZXMKICBzY2FsZV9hbHBoYV9jb250aW51b3VzKHJhbmdlID0gYyguMywgMSkpICsKICBzY2FsZV9zaXplX2NvbnRpbnVvdXMocmFuZ2UgPSBjKC41LCA4KSkgKwogICMgd2lkZW4gcGxvdHRpbmcgYXJlYQogIHhsaW0oTkEsIDEwMCkKYGBgCgpgYGB7cn0KIyBoaWdobGlnaHQgY2xhZGVzIHdpdGggYmFja2dyb3VuZCBjb2xvcnMKcDIgPC0gcCArIAogIGdlb21faGlsaWdodChub2RlID0gNDkzLCBmaWxsID0gY29sc1sxXSwgYWxwaGEgPSAuMykgKwogIGdlb21faGlsaWdodChub2RlID0gNDk2LCBmaWxsID0gY29sc1sxXSwgYWxwaGEgPSAuMykgKwogIGdlb21faGlsaWdodChub2RlID0gNDI5LCBmaWxsID0gY29sc1syXSwgYWxwaGEgPSAuMykgKwogIGdlb21faGlsaWdodChub2RlID0gMzAzLCBmaWxsID0gY29sc1szXSwgYWxwaGEgPSAuMykgKwogIGdlb21faGlsaWdodChub2RlID0gNDA4LCBmaWxsID0gY29sc1s0XSwgYWxwaGEgPSAuMykgKwogICMgcGxvdCB0cmVlIGFnYWluIHRvIGJlIG9uIHRvcCBvZiB0aGUgaGlnaGxpZ2h0cwogIGdlb21fdHJlZSgpICsKICBnZW9tX3Jvb3Rwb2ludChzaXplID0gMSkKYGBgCgpgYGB7cn0KIyBoaWdobGlnaHQgY2xhZGVzIHdpdGggYnJhbmNoIGNvbG9ycwpwMyA8LSBwICsgCiAgYWVzKGNvbCA9IGdyb3VwKSArCiAgc2NhbGVfY29sb3JfbWFudWFsKHZhbHVlcyA9IGMoImdyYXkzMCIsIGNvbHMpKQpgYGAKCmBgYHtyfQpwbG90cyA8LSBtZ2V0KGMoInAiLCAicDIiLCAicDMiKSkKYGBgCgpgYGB7ciwgZmlnLndpZHRoPTE4LCBmaWcuaGVpZ2h0PTZ9CmdyaWQuYXJyYW5nZShwLCBwMiwgcDMsIG5yb3cgPSAxKQpgYGAKCmBgYHtyfQojIHBuZyB3aXRoIDN4MQpnZ3NhdmUoIi4uL2dyYXBocy9waHlsb19mdWxsLnBuZyIsIGFycmFuZ2VHcm9iKGdyb2JzID0gcGxvdHMsIG5yb3cgPSAxKSwgCiAgICAgICB3aWR0aCA9IDI0LCBoZWlnaHQgPSA4LCBzY2FsZSA9IDIsIGRwaSA9IDcyKQoKIyBwZGYgd2l0aCAxIHBlciBwYWdlCmdnc2F2ZSgiLi4vZ3JhcGhzL3BoeWxvX2Z1bGwucGRmIiwgbWFycmFuZ2VHcm9iKGdyb2JzID0gcGxvdHMsIG5yb3cgPSAxLCBuY29sID0gMSksIAogICAgICAgd2lkdGggPSA4LCBoZWlnaHQgPSA4LCBzY2FsZSA9IDIsIGRwaSA9IDcyKQpgYGAKCiMgU2FtcGxlIHNpemUgaW4gZGV0YWlsCgpgYGB7cn0KIyBzdWJzZXQgdHJlZSB0byBqdXN0IHRob3NlIHNwZWNpZXMgd2hvIGhhdmUgc2FtcGxlIHNpemVzIHJlcG9ydGVkLCBpLmUuIHRob3NlIHdobyB3ZXJlIHRlc3RlZAp0b19kcm9wIDwtIHRyZWUzICU+JSBmaWx0ZXIoaXMubmEobikpICU+JSBwdWxsKGxhYmVsKQp0cmVlNSA8LSBkcm9wLnRpcCh0cmVlNCwgdG9fZHJvcCkKYGBgCgpgYGB7cn0KZGF0YTMgPC0gZGF0YSAlPiUgCiAgc2VsZWN0KGxhYmVsLCBldmVyeXRoaW5nKCkpICU+JSAKICByZW5hbWUobnVtID0gbikKCiMgc3BlY2llcyB3aXRoIG1vcmUgdGhhbiBYIHNpdGVzIGNhbiBnZXQgYSBkZW5zaXR5PwpkM2EgPC0gZGF0YTMgJT4lIGdyb3VwX2J5KHNwZWNpZXMpICU+JSBmaWx0ZXIobl9kaXN0aW5jdChzaXRlKSA+PSA0KQpkM2IgPC0gZGF0YTMgJT4lICMgc2V0ZGlmZihkYXRhMywgZDNhKSAjIyA8LSBkbyBzZXRkaWZmIGluc3RlYWQgdG8gTk9UIHNob3cgcG9pbnRzIGZvciBkZW5zaXRpZXMKICBncm91cF9ieShzcGVjaWVzKSAlPiUgCiAgIyBjcmVhdGUgdmFyaWFibGUgbnVtMiBpcyBOQSBpZiB0aGVyZSdzIG9ubHkgb25lIGRhdGEgcG9pbnQgZm9yIGEgc3BlY2llcwogICMgLS0+IHRob3NlIHNwZWNpZXMgd2lsbCBvbmx5IGdldCB0aGUgdmVydGljYWwgY3Jvc3NiYXIKICBtdXRhdGUoZmxhZyA9IG5fZGlzdGluY3Qoc2l0ZSkgPT0gMSkgJT4lIAogIHVuZ3JvdXAgJT4lIAogIG11dGF0ZShudW0yID0gaWZlbHNlKGZsYWcsIE5BLCBudW0pKQoKIyBmb3IgdmVydGljYWwgY3Jvc3NiYXIgPSBtZWRpYW4KZGF0YTQgPC0gZGF0YTMgJT4lIAogIGdyb3VwX2J5KGxhYmVsLCBzcGVjaWVzKSAlPiUgCiAgc3VtbWFyaXNlKE1kbiA9IG1lZGlhbihudW0pKSAjIHRvdGFsTiA9IHN1bShudW0pLCBzaXRlc04gPSBuX2Rpc3RpbmN0KHNpdGUpCmBgYAoKYGBge3J9CnEgPC0gZ2d0cmVlKHRyZWU1LCBhZXMoY29sID0gZ3JvdXApKSArCiAgIyB0aGlzIGlzIGEgZHVtbXkgcG9pbnQgdG8gZXhwYW5kIHRoZSB4IHNjYWxlCiAgIyB0aGUgdHlwaWNhbCB3YXlzIHdlaXJkbHkgYWxzbyBleHBhbmQgaXQgZm9yIHRoZSBzYW1wbGUgc2l6ZSBwYW5lbCBvbmNlIHRoYXQncyBhZGRlZAogIGdlb21fcG9pbnQoZGF0YSA9IHRpYmJsZSh4ID0gMTM1LCB5ID0gMSwgLnBhbmVsID0gIlRyZWUiLCBncm91cCA9IE5BKSkgKwogICMgdGlwIGxhYmVscwogIGdlb21fdGlwcG9pbnQoYWVzKHNpemUgPSBuKSwgYWxwaGEgPSAuNSkgKwogIGdlb21fdGlwbGFiKG9mZnNldCA9IDQsIHNpemUgPSAzKSArCiAgIyB0d2VhayBzY2FsZXMKICBzY2FsZV9jb2xvcl9tYW51YWwodmFsdWVzID0gYygiZ3JleTMwIiwgY29scykpICsKICBzY2FsZV9maWxsX21hbnVhbCh2YWx1ZXMgPSBjb2xzWzRdKSArICMgd2hlbiBhbGwgY2F0ZWdvcmllcyBhcmUgdGFrZW46IGNvbHMKICAjIGRpc3BsYXkgdGltZXNjYWxlIGF0IHRoZSBib3R0b20KICB0aGVtZV90cmVlMigpCmBgYAoKYGBge3IsIGZpZy53aWR0aD04LCBmaWcuaGVpZ2h0PTJ9CiMgdHJ5aW5nIG91dCB3aGljaCBnZW9tcyB0byBsYXllcgpzMSA8LSBnZ3Bsb3QoZGF0YTMsIGFlcyh4ID0gbnVtLCB5ID0gbGFiZWwsIGdyb3VwID0gbGFiZWwpKSArCiAgZ2VvbV9kZW5zaXR5X3JpZGdlcyhsd2QgPSAuMywgY29sID0gImdyZXk4MCIsIGJhbmR3aWR0aCA9IDEpCgpzMiA8LSBnZ3Bsb3QoZGF0YTMsIGFlcyh4ID0gbnVtLCB5ID0gbGFiZWwsIGdyb3VwID0gbGFiZWwpKSArCiAgZ2VvbV9jcm9zc2JhcmgoZGF0YSA9IGRhdGE0LCBhZXMoeCA9IE1kbiwgeG1pbiA9IE1kbiwgeG1heCA9IE1kbiksIHdpZHRoID0gLjUpICsKICBnZW9tX3BvaW50KGFscGhhID0gLjUpCgpzMyA8LSBnZ3Bsb3QoZGF0YTMsIGFlcyh4ID0gbnVtLCB5ID0gbGFiZWwsIGdyb3VwID0gbGFiZWwpKSArCiAgZ2VvbV9ib3hwbG90aCgpCgpwbG90X2dyaWQoczEsIHMyLCBzMywgbnJvdyA9IDEpCmBgYAoKYGBge3IsIGZpZy53aWR0aD00LCBmaWcuaGVpZ2h0PTN9CiMgcmlnaHQtc2lkZSB2aXogZGVwZW5kcyBvbiB0aGUgbnVtYmVyIG9mIHNpdGVzIHBlciBzcGVjaWVzOgojIDEgc2l0ZSA9IHZlcnRpY2FsIGNyb3NzYmFyIG9ubHkKIyAyKyBzaXRlcyA9IHBvaW50cyArIGNyb3NzYmFyIGF0IG1lZGlhbgojIFgrIHNpdGVzID0gZGVuc2l0aWVzIChjdXJyZW50bHksIFggPSA0IGp1c3QgdG8gaWxsdXN0cmF0ZSkKCiMgZGlydHkgaGFjazogdGhlcmUncyBhIHNtYWxsIGNoYXJhY3RlciBpbiBmcm9udCBvZiAiU2FtcGxlIHNpemVzIiB0byBoYXZlIHRoYXQgcGFuZWwgc29ydCB0byB0aGUgcmlnaHQgKGFscGhhYmV0aWNhbGx5KSB1bnRpbCBJIGZpZ3VyZSBvdXQgd2h5IGl0IGRvZXNuJ3QganVzdCBnbyBieSBvcmRlci4gVGhpcyBjcm9wcGVkIHVwIGFzIGFuIGlzc3VlIHdoZW4gSSBhZGRlZCB0aGUgZHVtbXkgcG9pbnQgZm9yIHRoZSB4LWF4aXMgZXhwYW5zaW9uLi4uCgpxICU+JSAKICBmYWNldF9wbG90KCLZtCBTYW1wbGUgc2l6ZXMiLCAKICAgICAgICAgICAgIGQzYSwgZ2VvbV9kZW5zaXR5X3JpZGdlcywgYWVzKHggPSBudW0sIGdyb3VwID0gbGFiZWwsIGZpbGwgPSBncm91cCksIAogICAgICAgICAgICAgYWxwaGEgPSAuNSwgbHdkID0gLjMpICU+JSAKICBmYWNldF9wbG90KCLZtCBTYW1wbGUgc2l6ZXMiLAogICAgICAgICAgICAgZGF0YTQsIGdlb21fY3Jvc3NiYXJoLCBhZXMoeCA9IE1kbiwgeG1pbiA9IE1kbiwgeG1heCA9IE1kbiwgZ3JvdXAgPSBsYWJlbCwgCiAgICAgICAgICAgICBjb2wgPSBncm91cCksIGFscGhhID0gLjUsIHdpZHRoID0gLjMsIHBvc2l0aW9uID0gcG9zaXRpb25fbnVkZ2UoeSA9IC0uMSkpICU+JQogIGZhY2V0X3Bsb3QoItm0IFNhbXBsZSBzaXplcyIsIAogICAgICAgICAgICAgZDNiLCBnZW9tX3BvaW50LCBhZXMoeCA9IG51bTIsIGdyb3VwID0gbGFiZWwpLCBhbHBoYSA9IC41LCAKICAgICAgICAgICAgIHBvc2l0aW9uID0gcG9zaXRpb25fbnVkZ2UoeSA9IC0uMSkpICsKICBwYW5lbF9ib3JkZXIoKQpgYGAKCmBgYHtyfQpnZ3NhdmUoIi4uL2dyYXBocy9waHlsb19yaWRnZS5wbmciLCB3aWR0aCA9IDQsIGhlaWdodCA9IDMsIHNjYWxlID0gMikKYGBgCgpSZWxldmFudCBmdW5jdGlvbnMgdG8gdGlua2VyIHdpdGgsIGlmIG5lY2Vzc2FyeToKCmBgYHtyfQpmYWNldF9wbG90CmBgYAoKYGBge3J9CmdndHJlZTo6OmFkZF9wYW5lbApgYGAKCiMgU2Vzc2lvbiBpbmZvCgpgYGB7cn0Kc2Vzc2lvbkluZm8oKQpgYGAKCg==